clear all
clc

cd('../../../../../cvx_2.2.2')
cvx_setup

warning('off','MATLAB:nargchk:deprecated')
warning('off','MATLAB:illConditionedMatrix')

cd('../bernstein_simulations/J2/codes/discrete_z')

%% DGP

J = 2;    % # of goods
N = 20000;  % # of individuals

alpha = 1;
beta = 1;

mu_x = 0;
sigma_x = 1;

nsupp_z = 10;
supp_z = linspace(0.1,0.9,nsupp_z);
increment_z = mean(diff(supp_z));

mu_c = -0.35;
sigma_c = 1.5;

% PREP FOR ESTIMATION

m_regr = [4 7 4*ones(1,J-2) 7 7 4*ones(1,J-2)];

[combos_all_unique] = combinations_2(m_regr,J);

% constraints

[Aineq,bineq] = constraint_mat_theta_2(combos_all_unique,J);

Aeq = [];
beq = [];

ntheta = size(Aineq,2);

ss=50;

rng(10101010*ss+ss^2)

tic

[data,data_trans] = generate_data_simult_categ(N,J,supp_z,alpha,beta,mu_c,sigma_c);

%% ESTIMATION ROUTINE

[theta_NPD,fval_NPD,flag_NPD] = estimation_NPD_cvx(data_trans,J,N,m_regr,Aineq,Aeq,bineq,beq,...
    combos_all_unique,ntheta);

toc

save('../../results/discrete_z/simult/MC_50.mat')

% generate data where z goes up to the penultimate point of support (to be
% able to increase values by one notch)
supp_z_reduced = supp_z(1:end-1);
rng(123)
[data,data_trans] = generate_data_simult_categ(N,J,supp_z_reduced,alpha,beta,mu_c,sigma_c);

script_ratio_categ

clear data* bineq Aineq input_list output_list

save('../../results/discrete_z/simult/MC_50.mat')




